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The WKB-approximation for the Bogoliubov-equations of 
the quasi-particle excitations in Bose-gases with condensate 
is worked out in the case of spherically symmetric trap po- 
tentials on the basis of the resulting quantization rule. The 
excitation spectrum is calculated numerically and also ana- 
lytically in certain limiting cases. It is found that the energy 
levels of a Bohr-Sommerfeld type quantization may be con- 
siderably shifted when the classical turning point gets close 
to the surface of the condensate. 
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The experimental realization and study of Bose- 
Einstein condensates in alkali atom gases confined by 
magnetic traps has induced a vivid activity in the 
theoretical investigation of such systems (See || for a 
recent review and for further references). From a the- 
oretical point of view the existence of the external po- 
tential requires new methods for calculating the physical 
properties of the quantum gases with Bose-condensation. 
Our aim is to solve the Bogoliubov-equations in WKB- 
approximation and to determine the excitation spectrum 
on the basis of the resulting quantization rule. 

In the Bogoliubov-theory the field operator can be ex- 
pressed as a linear combination of quasiparticle creation 
and annihilation operators. The corresponding (nonuni- 
form) expansion coefficients Uj(r) and Vj(r) obey the 
coupled linear Bogoliubov eigenvalue equations JtJ 



Hhf 
-K*(r) 



-K(r) 
Hhf 



-Vj(r) 



(1) 



where j denotes one of the quasiparticle states and Ej 
is the corresponding quasiparticle energy. The Hartree- 
Fock operator Hhf takes the form 



H 



HF 



- — A + U(r) + 2\K(r)\ 
2m 



(2) 



where U(r) is the trap potential, /i is the chemical po- 
tential, 



K(r) = 



Mr) 2 



(3) 



denotes the potential-like contribution of the con- 
densate, whose wave- function ipo( r ) is normalized as 



/ d 3 r |-0o( r )| 2 = Nq. Nq is the number of particles in 
the condensate and a is the s-wave scattering length. In 
the following we shall assume that a > 0. The quasiparti- 
cle amplitudes Uj(r) and Vj(r) are normalized according 
tofl 



d 3 r (u*(r)u k (r) - v*{r)v k {r)) = S jk 



(4) 



For the sake of simplicity we choose the external potential 
as spherically symmetric. Moreover we shall take ipo{r) 
and hence also K(r) as real and shall also make frequent 
use of the Thomas-Fermi approximation |^] , which leads 
to 



|^o(r)| 2 = (4^(^-^W) if r<r T F 
1 otherwise 



(5) 



Here, U{ttf) = Mtf and P"TF is fixed by normalization. 

One can introduce spherical coordinates r, 9, <f> and 
separate variables in the usual way: 



(r)\ _ 1 (u nl (r) 



Vnl(r) 



Yi, 



(6) 



where j denotes the usual quantum numbers (n, I, m) for 
isotropic problems and the Yi m are the spherical harmon- 
ics. 

To solve the coupled, radial equations obtained from 
(|l|) it is advantageous to use the linear combinations || 

G nl( r ) = {Unl{r)±V n l{r)) , (7) 

which satisfy the uncoupled equations 

' H HF - K(r) 2 -E 2 T [H HF ,K(r)}) G ± (r) = . (8) 



Here [, ] denotes the commutator. (For brevity we have 
omitted the indices n and I). Furthermore, it follows 
from the original equations that 



G ± = ]-[H H F±K{r))GT 



(9) 



which is compatible with equations (||). 
Now the operator Hhf has the form 



H 



HF 



2m dr 2 



+ U ef f(r), 



(10) 



1 



where 

U eff(r) = + U W + 2K(r) (11) 

In our WKB treatment we use Langer's rule by replacing 
1(1+1) by (1+1/2) 2 . In the following U e ff (r) is considered 
as a classical potential. 

We shall consider two types of solutions 



G+( r ) =e *( So+ * Sl+ -), 
G+(r) = e -^(So+SS 1+ ...) 



(12) 
(13) 



with real functions S'o(r), S'i(r), . . . and So(r), Si(r), . . . 
respectively. Gathering terms having different powers of 
h one gets first order ordinary differental eqations for the 
unknown quantities occuring in (|l2] ) and (O). 

First we consider solutions of the form (|l2|). The 0(h) 
equation is the classical Hamilton- Jacobi equation for the 
radial action So (r) , from which one can express the clas- 
sical radial momenta as 



\Pr\ 



dS Q 



dr 



2m(±y / E 2 + K 2 -U, 



eff 



(14) 



We shall assume that U e ft > in which case only 
the plus sign is allowed to have p r real. This is the 
case for instance in the Thomas-Fermi approximation 
(|J). We introduce the radial velocity in the usual way 
v r = dH/dp r by regarding E in (Jl4|) as the classical 
Hamiltonian H(p r ,r). The obtained expression JIc| ] 



^/E 2 + K 2 (r) Pr 



E 



(15) 



reflects the peculiarity of the classical quasi particle dy- 
namics in traps. The effective quasi-particle mass, which 
can be read off from (|l4|), is energy and space-dependent. 
It approaches the particle mass at the boundary of the 
condensate, but can become much smaller yet remains 
non-zero even in the center of very large condensates in 
traps. This is a fundamental difference to the untrapped 
case, where the limit E — > can be taken, in which the 
quasi-particle mass vanishes. 

By solving (|8|) with the ansatz ( |l2| ) up to Si, then 
using (Q) for G~(r), and finally transforming back from 
G ± (r) to u n i(r) and v n i(r), particular solutions of the 
radial Bogoliubov equations are obtained in the form 



u„i(r) 
v„i (r) 



where u 2 n = 



Const x 



u B (r) 
v B (r) 



1 



,±£ J r Pr(r) dr 



(y/l + (K(r)/E)* + l) /2, 



(16) 



1 are 



the generalizations of the usual Bogoliubov-coefficients 
for the case without trapping potential. Note that the 
classical probability distribution is inversely proportional 
to | iv | as expected physically. 



Solutions ( fig ) are valid in the classically allowed re- 
gion, i.e., between the classical turning points r t \ and 
Tt2 > f t i defined by the condition p r (ru) = 0, i = 1,2. 
We shall assume that there are two turning points only. 
There may be three cases. Case A: if r t \ < ttf < Tti, in 
other words, the classical particle enters the condensate, 
then leaves it, and returns back again etc.; Case B: if 
ttf < i~ti < r t 2, i.e., we have only a simple classical mo- 
tion in the trapping potential; Case C: if r t \ < r t 2 < rxF, 
in which case the classical motion is confined to the con- 
densate. 

Next we construct solutions of the form (O) proceed- 
ing similarly as before. Using the ansatz (|l3|) in Eq. (||). 
there can exist two different solutions for So, 



,(<) 



dS t 



dr 



1,2. 



(17) 



Both signs are allowed for example in the Thomas-Fermi 
approximation (^). The solution for i = 1 is defined only 
outside the classically accessible region, while the other 
one (i = 2) is permissible for all r- values, if U e ff > (as 
we suppose), and represents a solution which can only 
occur in the two component quasi particle dynamics. Let 
us define furthermore quantities wi^ similar as in (|l5| ) by 
the relations 



■,W = 



v/£ 2 + K 2 (r) q® 



E 



(18) 



For a smooth potential U e f / (r) it can be shown that nor- 

~ (2) 

malizable eigenfunctions cannot contain S in a WKB 
solution. However for states whose radial wavelength 
near the characteristic radius of the condensate is large 
compared to the width of the boundary layer there, the 
effective potential can no longer be treated as smooth. 
Furthermore, such a type of contribution must always 
be present in v(r) asymptotically, if the condensate is 
restricted to a finite region. 

Let us now consider the allowed solutions in case A. 
Requiring normalizibily and performing turning point 
matching at r t j one obtains 



Unl 

v n i 





( v B 




\-u B 


vW 2) i 






(u B \ 



V B 



exp 



F(r), 



qW(r) dr 
(19) 



where j = 1 and j = 2 correspond to r < ttf and 
r > ttf respectively. C\j and are arbitrary constants 
and 



F(r) = exp 



Z r = w { }\ 



(-iy 



for 



q { r l) {r) dr 

< r < r t i 
r t 2 < r 



(j 

u 



1) , 

2) , 



(20) 



2 



and 

F(r) = 2 sin 



i-iy 



h 

for 



p r (r) <ir 



r tl <r < r TF (j = 1), 
r TF <r <r t2 (j = 2). 



(21) 



Requiring that u(r), v(r) and their first derivatives 
are continuous at ttf one gets four homogeneous lin- 
ear equations for the four unknown constants. In order 
to get non-trival solutions the determinant of the coeffi- 
cient matrix should vanish. This leads to the scmiclassi- 
cal quantization rule 



= — — — cos 
n 



mL 



Ia 



Pa 



h 

L 
2E 



' Sin (l 



TT\ . ( Ib , 7T 

4 ri + 4 



hmL 



(22) 



where we have introduced the notations pa = Pr(fTF), 
q { r 2) (r TF ), L = {dK/dr) rTP+0 - (dK/dr) r 



PB 



)r T F—0> 



l A = f r ™ Pr dr, I B = J r " F p r dr. 

Keeping only the first term on the right hand side of 
(p2|) leads to the usual Bohr-Sommerfeld quantization 
rule 



(23) 



with the integer radial quantum number n > and in- 
cluding the Maslov-indices due to the two turning points 
in the radial motion. The case B can be treated in an 
analogous way leading to the quantization rule (^3j) with 
K{r) = within the range of integration. 

We discuss first the energy levels on the basis of ( |23| ) 
and will turn back to the consequences of the general ex- 
pression ( p2|) afterwards. To evaluate (23) we choose a 
harmonic potential U(r) = mco^r 2 /2 often used in theo- 
retical consideration pj^,p|,pl| . To distinguish between 
cases A and B let us use the dimensionless variables 
J = huJo(l + l/2)/(2/z), E — E/fi. For energies and 
angular momenta in the region 1<2J — 1 < E < J 2 
case B occurs and the energy spectrum is simply that of 
a harmonic oscillator shifted by /x 



E 



(osc) 
n,l 



Tiluq [ 2n + / + 



(24) 



The self-consistency condition for case B is then 1 + 1/2 > 
2/j,/Hujo + \jA:[i(2n + 1)/Twjq . One can check that case C, 
i.e., the classical motion is entirely inside the condensate, 
is not possible. 

Considering the nontrivial case A in region < J 2 < E 
the action integral in Eq. (|23|) can still be performed 



analytically (For the details see Ref. |i2[|), but the re- 
sult is rather cumbersome, and the energies E n i can- 
not be expressed explicitely. However, due to (|23|) the 
semiclassical energies fulfill the scaling relation E n j = 
hu>oG ni i(fJ./hujo)- We discuss here some limiting cases. 
One interesting limit is when one considers the high-lying 
levels, i.e., when E n j 3> /i is fulfilled. Then the main con- 



tribution to the action integral in (23) comes from the 



region outside of the condensate, leading to a spectrum 
which is almost that of a harmonic oscillator. Expanding 
the action integral to the next to leading correction in 
(Jt/En i one gets: 



E n ,i = £ , i° SC) + ^>o5 n ,i 



3/2 



3tt 



(25) 



This result and that of the perturbation theoretical cal- 
culation JlU agree for large (2n + I) values. 

The other interesting limit is the region of excitation 
energies small compared with the chemical potential. To 
reach it formally, the angular momentum and the radial 
quantum numbers / and n are kept fixed but /i/huo tends 
to infinity. The main contribution to the radial action 
integral ( p3|) comes from those r values, which are within 
the condensate. To leading order: 



E U: i ~ huj Q [2n 2 + 2nl + 'in + I + l] 



1/2 



(26) 



Our result (|26|) almost coincides with that of Stringari's 
hydrodynamic calculatio n |llf| , except for the last con- 
stant 1 within [...J 1 / 2 in (|26|), which has an appreciable 
effect only on the lowest levels. It is, of course, not unex- 
pected that a WKB approach may fail there. For some- 
what higher energies at fixed but large chemical poten- 
tial there is a considerable region where the two spectra 
calculated in WKB and in hydrodynamical approxima- 
tions, respectively, overlap. For even higher energies the 
applicability of the hydrodynamical approach looses its 
validity. The task of solving ( p3| ) for E can be carried out 
numerically in a straightforward manner for given scat- 
tering length a, trapping potential (i.e. ojq) and number 
A^o of atoms in the condensate, fixing the single param- 
eter 2h/TujJo = (15N a/ \J%/muj ) 2 / b on which the spec- 
trum depends. An example of the results obtained is 
depicted in Figure |l|. 

Let us turn now to the discussion of (p2|). By solv- 



ing (22) numerically for the energies we have found that 
the corrections AE n j to the levels defined by (|||) are 
small except when the classical inner turning point r t \ 
gets close to the surface of the condensate (the border 
between regions A, B in the (n,l)-plane), in which case 
the radial wavenumber pa goes to zero. A.E n j then be- 
comes large but decreases rapidly when going away from 



3 



this situation. Note that in this case the classical or- 
bits are just glancing at the surface of the condensate. 
The solution of the Bogoliubov equations then contains 
an anomalous contribution, namely, the first term on the 
right hand side of (|l9|), which is exponentially localized 
at the surface. The effect remains even for high ener- 
gies when the anomalous part of u and the usual part 
of v become negligible. At such energies it is generally 
assumed that the Bogoliubov equations go over to the 
Hartree-Fock equations. Our results suggest that there 
are exceptional states at the border of region A and B for 
which this is not true due to the disturbance at the sur- 
face. Physicallythe effect is caused by the narrow bound- 
ary layer of the condensate which looks effectivly sharp 
for orbits glancing on the surface. Its qualitative aspects 
can therefore be expected to be independent from the 
WKB and Thomas-Fermi approximations. Experimen- 
tally such anomalous states could be observed by their 
excitation via modulations of the trapping potentials as 
m @L @ or by light scattering. 

In this paper we have restricted ourselves, for the sake 
of simplicity, to the case of the spherically symmetric trap 
potential. Calculations along these lines for anisotropic 
harmonic oscillator trap potentials as they are used in the 
experiments will be published in a separate paper. 
Here we only mention that the corresponding classical 
Hamiltonian shows chaotic behavior Q], especially for 
energies comparable to the chemical potential. 
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FIG. 1. Energy levels E(n,l) = E n i/hu>o obtained numer- 
ically from ( p3| ) for an isotropic harmonic oscillator trap po- 
tential U(r) = mujQr 2 /2 as a function of the radial quantum 
number n and of the angular momentum quantum number I. 
The chemical potential was chosen to be 8huo 
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